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Ab^ract 

In  this  paper,  we  present  a  method  for  recovering  both  the  shape  of  an  object  and  its  motion  relative  to 
the  camera  from  a  sequence  of  images  of  the  object,  using  feature  points  tracked  throughout  the 
sequence.  The  method  is  based  on  the  factorizmion  method  developed  by  Tomasi  &  Kanade,  which 
achieves  its  accuracy  by  using  a  large  number  of  pdnts  and  images,  and  robustly  applying  a  well-under¬ 
stood  matrix  computation,  the  singular  value  decomposition,  to  the  highly  redundant  input  data. 

While  the  Tomasi  &  Kanade  method  was  based  on  orthographic  projection,  our  method  uses  a  projec¬ 
tion  model  known  as  parapospective  projection.  Orthographic  projection  does  not  account  for  the 
apparent  change  in  size  of  an  object  as  it  moves  toward  or  away  from  the  camera,  nor  the  different  angle 
from  which  an  object  is  viewed  as  it  moves  parallel  to  the  image  plane.  In  contrast,  paraperspective  pro¬ 
jection  closely  ^jptoximates  perspective  projection  by  modelling  both  of  these  effects.  A  new  formula- 
tirm  based  on  this  projection  mo^I  allows  us  to  apply  the  bctorization  method  to  a  wider  range  of 
scenarios,  and  to  recover  the  distance  from  the  camera  to  the  object.  The  method  assumes  no  model  of 
the  motion  or  of  the  object’s  shape,  arxl  recovers  the  shape  and  motion  accurately  even  fcx  distant 
objects. 

We  present  several  experiments  which  illustrate  the  method’s  performance  ova-  a  wide  range  of  noise 
values  and  a  wide  range  of  distances  from  the  object 
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1.  Introduction 

Recovering  the  geometry  of  a  scene  and  the  motion  of  the  camera  from  a  stream  of  images  is 
an  important  task  in  a  variety  of  applications,  including  navigation,  robotic  manipulation, 
and  aerial  cartography.  While  this  is  possible  in  principle,  traditional  methods  have  failed  to 
produce  reliable  results  in  many  situations  [2]. 

Tomasi  and  Kanade  [7][8]  developed  a  robust  and  efficient  method  for  accurately  recovering 
the  shi^  and  motion  of  an  object  from  a  sequence  of  images  using  extracted  feature  points. 
Their  method  uses  an  orthographic  projection  model,  which  is  described  by  linear  equations. 
It  achieves  its  accuracy  and  robusmess  by  using  a  large  number  of  frames  and  feature  points, 
and  by  directly  computing  shape  without  computing  the  depth  as  an  intermediate  step.  The 
method  was  tested  on  a  variety  of  real  and  synthetic  images,  and  was  shown  to  perform  well 
even  for  distant  objects. 

There  are,  however,  some  limitations  of  the  method  due  to  its  use  of  the  orthographic  pro¬ 
jection  model.  The  model  contains  no  notion  at  all  of  the  distance  from  the  camera  to  the 
object  As  a  result,  image  sequences  containing  large  translations  toward  or  away  from  the 
camera  often  produce  deformed  object  shapes,  as  the  method  tries  to  explain  the  size  differ¬ 
ences  in  the  images  by  creating  size  differences  in  the  object  It  also  supplies  no  estimation 
of  translation  along  the  camera’s  optical  axis,  which  limits  its  usefulness  for  certain  tasks. 

Fortunately,  there  exist  several  perspective  approximations  which  capture  more  of  the 
effects  of  perspective  projection  while  remaining  linear.  Scaled  orthographic  projection, 
sometimes  referred  to  as  “weak  perspective”  [3],  accounts  for  the  scaling  effect  of  an  object 
as  it  moves  towards  and  away  from  the  camera.  Paraperspective  projection,  first  introduced 
by  Ohta  in  [4]  and  named  by  Aloimonos  in  [1],  accounts  for  the  scaling  effect  as  well  as  the 
different  angle  from  which  an  object  is  viewed  as  it  moves  in  a  direction  parallel  to  the 
image  plane. 

In  this  paper,  we  present  a  new  factorization  method  based  on  the  paraperspective  projection 
model.  The  paraperspective  factorization  method  takes  a  set  of  points  extracted  from  the 
images  and  tracked  from  one  image  to  the  next,  and  computes  the  shape  of  the  object  and 
the  motion  of  the  camera.  The  method  assumes  no  model  of  the  motion  or  of  the  object’s 
shape,  and  is  still  fast  and  robust  with  respect  to  noise.  However,  it  can  be  applied  to  a  wider 
realm  of  situations  than  the  original  factorization  method,  such  as  sequences  containing  sig¬ 
nificant  depth  translation  or  containing  objects  close  to  the  camera,  and  can  be  used  in  appli¬ 
cations  where  it  is  important  to  recover  the  distance  to  the  object,  such  as  navigation. 

We  begin  by  describing  our  camera  and  world  reference  frames  and  introduce  the  mathe¬ 
matical  notation  that  we  use.  We  review  the  original  factorization  method  as  defined  in  [8], 
presenting  it  in  a  slightly  different  manner  in  order  to  make  its  relation  to  the  paraperspec¬ 
tive  method  more  apparent  We  then  present  our  paraperspective  factorization  method.  We 
conclude  with  the  results  of  some  experiments  which  demonstrate  the  practicality  of  our 
system. 
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2.  The  Problem 

In  a  shape-£rom-motion  pioblens,  we  ate  given  a  sequence  of  F  images  taken  from  a  camera 
that  is  moving  relative  to  an  object  Imagine  that  we  locate  P  prominent  feature  points  in  the 
first  image,  and  track  these  points  from  each  ima^  to  the  next  recording  the  coordinates 
of  each  point  p  in  each  image  f.  Each  feature  point  p  that  we  track  corresponds  m  a 
single  world  point  located  at  position  in  some  fixed  world  coordinate  system.  Each 
image  f  was  t^en  at  some  specific  camera  orientation,  which  we  describe  by  the  orthonor¬ 
mal  unit  vectors  ip  jp  and  where  points  along  the  camera’s  line  of  sight  y  corre¬ 
sponds  to  the  camera  image  plane’s  x-axis,  and  corresponds  to  the  camera  image’s  y-axis. 
There  are  of  course  only  three  independent  parameters  -  roll,  pitch,  and  yaw  -  to  describe  the 
camera’s  orientation,  and  the  camera  frame  vectors  can  be  written  in  terms  of  these  three 
variables.  Additionally,  we  describe  the  position  of  die  camera  in  each  frame  /  by  the  vector 
tf  pointing  to  the  camera’s  focal  point  Ihis  formulation  is  illustrated  in  Figure  1. 


WOTld 

Origin 


length/ 

F^nre  1 

Coordinate  System 


The  result  of  the  feature  tracker  is  set  of  P  feature  point  coordinates  for  each  of  the  F 

frames  of  the  image  sequence.  From  this  information,  our  goal  is  to  recover  the  estimated 
shape  of  the  object,  given  by  the  position  of  every  point,  and  the  estimated  motion  of  the 
camera,  given  by  i/,  j/,  k/  and  t/  for  each  frame  in  the  sequence.  Rather  than  recover  tf  in 
world  coordinates,  we  generally  recover  the  three  separate  components  t/  i/,  t/  j/,  and 

if  kf. 
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3.  The  Orthographic  Factorization  Method 

This  section  presents  a  sunomaiy  of  the  orthographic  factorization  method,  which  was  devel¬ 
oped  by  Tomasi  and  Kanade  in  1990.  A  more  detailed  description  of  the  method  can  be 
found  in  [7]  or  [8]. 

3.1.  Orthographic  Projection 

The  orthographic  projection  model  assumes  that  all  rays  are  projected  from  the  object  point 
parallel  to  the  camera’s  line  of  sight  so  that  they  strike  the  image  plane  orthogonally,  as 
illustrated  in  Figure  2. 

Image 


Figiire2 

Orthographic  Projection  in  two  dimensions 
Dotted  lines  indicate  true  perspective  projection 


Under  orthographic  projection,  a  point  p  whose  location  is  s^  will  be  observed  in  frame  f  at 
image  coordinates 

We  can  rewrite  these  equations  as 

Xfp  =  »/•  +  exf  yf^  =  »/■  +  cyr  (2) 

where 
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3.2.  Decomposition 

We  organize  all  of  the  feature  point  coordinates  into  &2FxP  measurement  matrix  w 

such  that 
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Each  column  of  the  measurement  matrix  contains  all  the  observations  for  a  single  point, 
while  each  row  contains  all  the  observed  x-cooidinates  or  y-coordinates  for  a  single  frame. 
We  can  combine  equation  (2)  for  all  points  and  all  frames  into  the  single  matrix  equation 


iy  =  MS  +  r[i  ...  i].  (6) 

where  M  is  the  2Fx3  motion  matrix,  5  is  the  3x/»  shape  matrix,  and  r  is  a  2Fx  i  transla¬ 
tion  vector. 


Up  to  this  point  we  have  not  put  any  restrictions  on  the  location  of  the  world  origin,  except 
that  it  be  stationary  with  respect  to  the  object  For  simplicity,  we  set  the  world  origin  at  the 
center-of-mass  of  the  object  denoted  by  c,  so  that 


c  = 


\ 

P 


(7) 


This  enables  us  to  compute  the  i'*  element  of  the  translation  vector  T  directly  from  W,  sim¬ 
ply  as  the  average  of  the  row  of  the  measurement  matrix.  We  then  subtract  out  the  trans¬ 
lation  from  W,  leaving  us  with  a  “registered”  matrix  W’.  Because  W'  is  the  product  of  a 
2Fx3  motion  matrix  M  and  the  3  x  p  shape  matrix  5,  it’s  rank  is  at  most  3.  We  use  singular 
value  decomposidon  to  factor  W  into 


IT  =  IkS.  (8) 

where  the  product  /kS  is  the  best  possible  rank  three  approximation  to  IT,  in  the  sense  of 
minimizing  the  sum  of  squares  difference  between  corresponding  elements  of  W  and 
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33,  Normalization 

The  decomposition  of  equation  (8)  is  not  unique.  In  fact,  any  3x3  non-singular  matrix  A 
and  its  inverse  could  be  inserted  between  A/  and  S,  and  their  product  would  still  equal  W . 
Thus  the  actual  motion  and  shape  are  given  by 


S  =  A''S' 


(9) 


when  the  3  X  3  invertible  matrix  A  is  selected  appropriately.  We  determine  this  matrix  A  by 
observing  that  the  motion  matrix  M  must  be  of  a  certain  form.  Because  y  and  are  unit  vec¬ 
tors,  we  derive  from  equation  (4)  that 


l"/->  I*/-'- 

and  because  they  are  orthogonal. 


11^=0.  (11) 

Equations  (10)  and  (11)  give  us  3F  equations  which  we  call  the  metric  constraints.  These 
constraints  have  the  property  that  they  are  linear  in  the  6  unique  elements  of  the  symmetric 
positive  definite  3x3  matrix  Q  =  A^A,  making  the  problem  a  simple  linear  data-fitting  prob¬ 
lem  which  can  be  easily  solved  for  its  least  sum-of-squares  error  solution.  Once  the  symmet¬ 
ric  matrix  Q  has  been  determined,  as  long  as  it  is  a  positive  definite  matrix,  it  is  a  simple 
matter  to  find  A. 

3.4.  Shape  and  Motion  Recovery 

Once  the  matrix  A  has  been  found,  the  shape  and  motion  can  easily  be  computed  from  equa¬ 
tion  (9).  In  the  case  that  the  resulting  motion  matrix  M  still  does  not  exactly  satisfy  the  met¬ 
ric  constraints  (since  we  have  only  solved  for  the  A  which  provides  the  best  fit  to  the  metric 
constraints),  we  simply  normalize  each  row  of  At  and  compute  the  orthonormal  1/  and  ]/ 
which  are  closest  to  our  computed  values.  There  is  an  ambiguity  in  the  solution  due  to  an 
arbitrary  world  coordinate  orientation.  We  eliminate  it  by  ali^fing  the  world  axes  with  the 
first  frame’s  camera  axes,  so  that  li  =  [i  o  and  ji  -  [o  i  ^  . 

3.5.  Performance 

The  orthographic  factorization  method  has  been  tested  on  a  variety  of  real  image  sequences 
and  has  been  shown  to  perform  extremely  well.  In  a  laboratory  experiment  in  which  the 
carr^ra  motion  was  carefully  controlled  and  recorded,  the  method  recovered  the  rotation  of 
the  camera  with  a  maximum  error  of  less  than  0.4  degrees,  and  with  an  average  error  far 
less.  The  object  shape  was  also  recovered  accurately,  with  errors  between  the  measured  and 
computed  distances  between  pairs  of  points  generally  less  than  1  percent 

The  method  was  also  tested  on  an  outdoor  scene,  in  which  a  hand-held  camera  was  used  to 
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view  the  comer  of  a  house.  Despite  the  obvious  jerky,  non-smooth  motion  of  the  camera,  the 
method  was  able  to  recover  the  shape  of  the  house  and  the  motion  of  the  camera  [8]. 
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4.  Paraperspective  Factorization  Method 


4.1.  Paraperspective  Projection 


Under  orthographic  projection,  as  an  object  translates  along  the  camera’s  optical  axis,  its 
size  in  the  image  remains  constant,  and  as  an  object  translates  parallel  to  the  image  plane,  its 
image  simply  translates  as  well  Under  perspective  projection,  however,  an  object’s  size 
changes  as  it  translates  along  the  camera’s  optical  axis;  furthermore,  as  an  object  translates 
parallel  to  the  image  plane,  the  image  both  translates  and  presents  a  slightly  Afferent  view 
of  the  object,  making  the  object  appear  to  rotate.  The  former  effect  is  the  “scaling  effect’’, 
while  the  latter  is  sometimes  referred  to  as  the  “position  effect”  of  perspective  projection 
[1],  because  the  amount  of  apparent  rotation  depends  on  the  object’s  position  in  the  image 
relative  to  the  center  of  projection. 

In  this  paper,  we  use  an  approximation  to  perspective  projection  known  as  paraperspective 
projection,  which  was  int^uced  by  Ohta  in  order  to  solve  a  shape  from  texture  problem. 
Paraperspective  projection  closely  approximates  perspective  projection  by  modelling  both 
the  scaling  effect  and  the  position  effect,  while  retaining  some  of  the  linear  properties  of 
orthographic  projection.  The  paraperspective  projection  of  an  object  onto  an  image,  illus¬ 
trated  in  Figure  3,  involves  two  steps. 

1.  The  points  of  the  object  are  projected  along  the  direction  of  the  line  connecting  the  focal  point  of 
the  camera  to  the  object’s  centa-of-mass,  onto  a  plane  parallel  to  the  image  plane  and  passing 
through  the  object’s  center-of-mass. 

2.  These  points  are  then  {Rejected  onto  the  ima^  plane  using  {xn:s{)ective  {>rojection.  Because  the 
points  are  all  on  a  plane  parallel  to  the  image  pltuie,  this  is  equivalent  to  simply  scaling  the  image 
by  the  ratio  of  the  camera  focal  length  and  the  distance  along  the  camera’s  optical  axis  to  the 
object’s  center-of-mass.^ 


In  general,  the  projection  of  a  point  p  along  direction  r,  onto  the  plane  with  normal  n  and 
distance  from  the  origin  d,  is  given  by  the  equation  p’  =  p  -  r.  Figure  4  illustrates  the 

first  step  of  paraperspective  projection.  We  project  the  object  point  s  along  the  direction 
c  -  tf,  which  is  the  direction  from  the  camera’s  focal  point  to  the  object’s  center-of-mass, 
onto  the  plane  defined  by  normal  and  distance  from  the  origin  c  •  giving 


(s^k^)-(ck.) 

~  (c-t^)  k^  (c-y- 


(12) 


The  perspective  projection  of  these  points  onto  the  image  plane  is  given  by  subtracting 
from  s'fp  to  give  the  position  of  the  point  in  the  camera’s  coerdinate  system,  and  then  scaling 
the  result  by  the  ratio  of  the  camera’s  focal  length  to  zp  the  depth  to  the  object’s  center-of- 
mass,  where  (c-y  kp  This  yields  the  coordinates  of  the  projection  in  the  image 


i.  The  fcaied  aMbognphic  projectioa  model  (atn  known  at  'wedc  penpective'^  ii  limilv  to  panpenpective  projection, 
except  that  the  diieclion  of  the  initial  projectioa  if  parallel  to  the  camera’s  optical  axis  rather  than  parallel  to  the  line  connecting 
the  object's  cemer-ttf-mass  to  the  camera's  focal  poim.  This  model  capuues  the  scaling  effea  of  perspective  projection,  but  not 
the  positian  effea.  For  disctttsioa  of  scaled  onhagnphk  projection,  tee  Appendix  L 
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plane. 


*/i»  “  V^* 

Substituting  (12)  into  (13)  and  simplifying  yields  the  general  paraperspective  equations  for 

Xfp^^yfp 


I  r  L-  (c-U  1 

^/P  =  Z^  I  •/-  z.  M  •  (^P  -«)  +  («-  V>  •  V> 

^  ^  (14) 

/  r  1 

yfp=7/\jr  ^^  — k^J-(s^-c)  +  (c-t^)  .j^} 

For  simplicity,  we  assume  unit  focal  length,  /  =  1 . 

We  have  not  up  to  this  point  put  any  requiren»nts  on  our  world  coordinate  system  except 
that  it  be  stationary  with  respect  to  the  object.  Thus  we  see  that  we  can  simplify  our  equa¬ 
tions  by  placing  the  world  origin  at  the  object’s  center-of-mass,  or  setting  c  =  0,  thus  reduc¬ 
ing  (14) to 


»,-<•/•  J/>> 

These  equations  can  be  rewritten  as 

where 


(15) 


(16) 


(17) 


ex.  *  - 


vv 


cy.  =  —!—L 

Zf  f  Zf 


'/ - 77~- 


(18) 

(19) 


Notice  that  equation  (16)  is  identical  to  its  counterpart  for  orthographic  projection,  equation 
(2),  although  the  corresponding  definitions  of  cx^  cy^  ny,  and  differ  somewhat.  This  sim¬ 
ilarity  enables  us  to  perform  the  basic  decomposition  of  the  matrix  in  exactly  the  same  man¬ 
ner  as  we  did  for  orthographic  projection. 


4.2.  Decomposition 

We  can  combine  equation  (16),  for  all  points  p  firom  l  to  P,  and  all  frames  f  from  l  to  F, 
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into  the  single  matrix  equation 


[l  ...  l],  (20) 


or 

!♦'  =  A/5  +  r[l  ...  l].  (21) 

where  ly  is  the  2F  x  p  measurement  matrix,  M  is  the  2F  x  3  motion  matrix,  S  is  the  3  x  F 
shape  matrix,  and  T  is  the  2F  x  i  translation  vector.  We  have  set  c  =  0,  so  by  definition 


■*11  •••  *l/» 

“l 

CXi 

. 

... 

... 

*F1  •••  *FF 
I'll  -  yip 

3S 

“f 

“l 

[*1  — 

CXp 

^yi 

ypi  •••  ypp 

."f 

cyp 

c  =  ^  I  Sp  =  0- 

P-i 


Using  this  and  equation  (16)  we  can  write 


p  p  p 

I  ^/p=  I  (m/  Sp  +  cx,)  =111,.  X  Sp  +  Pcx^=PcXf 

p  p  P 

S  >/p=  I  O/  Sp  +  cy,)  =n,.  X  ip  +  Pcyf^Pcyf 

psi  p*i 


Therefore  we  can  compute  cx,  and  cy,  immediately  from  the  image  data  as 


(22) 


(23) 


1  1  ^ 

^yf=pl.yfp- 

p " 1  p* 1 


(24) 


We  then  subtract  these  values  from  the  corresponding  rows  in  W,  giving  the  registered  mea¬ 
surement  matrix 


1*11  • 

•  *1F 

C*1 

®i 

*F1  • 

>11  • 

•  *FP 

•  yip 

- 

CXp 

cy, 

[i ...  i] . 

nip 

“l 

ypi  • 

■  ypp. 

cyp 

."f 

(25) 


Since  vr  is  the  product  of  two  matrices  each  of  rank  at  most  3,  W'  has  rank  at  most  3,  just  as 
it  did  in  the  orthographic  projection  case.  If  there  is  noise  present,  the  rank  of  IT  will  not  be 
exactly  3,  but  by  computing  the  Singular  Value  Decomposition  (SVD)  of  w  and  only 
retaining  the  largest  3  singular  values,  we  can  factor  W  into 
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w  =  JhS,  (26) 

where  A/  is  a  2F  X  3  matrix  and  S  is  a  3  x  />  matrix.  Using  the  SVD  to  perform  this  factoriza- 
ticn  guarantees  that  the  product  is  the  best  possible  rank  3  approximation  to  W ,  in  the 
sense  that  it  minimizes  the  sum  of  squares  difference  between  corresponding  elements  of  W' 
and  AVI. 


4  Normalization 


Just  as  in  the  orthographic  case,  the  decomposition  of  W  into  the  product  of  M  and  ^  is  not 
unique.  Any  3x3  matrix  A  and  its  inverse  could  be  inserted  between  AV  and  5,  and  the  prod¬ 
uct  would  still  equal  W .  We  need  to  find  the  matrix  A  that  gives  the  true  shape  and  motion 


=  kA 
=  A-‘S 


Again,  we  determine  this  matrix  A  by  observing  that  the  motion  matrix  M  must  be  of  a  cer¬ 
tain  form.  As  in  the  orthographic  case,  we  take  advantage  of  the  fact  that  and  are  unit 
vectors  and  are  orthogonal.  According  to  equation  (19),  we  observe  that 


(28) 


We  know  the  values  of  cx^  and  cjy  from  our  initial  registration  step,  but  because  we  do  not 
know  the  value  of  the  depth  zp  we  cannot  impose  individual  constraints  on  the  magnitudes 
of  and  as  we  did  in  the  orthographic  factorization  method.  However,  from  equation 
(28)  we  see  that 


I  .  I"/  _  !■/ 

1  +  1  +  cy^ 

Therefore  we  adopt  the  following  constraint  on  the  magnitudes  of  and  n^: 


(29) 


I"/  I"/  _  0 

1  +  cx^  1  +  cy^ 


(30) 


In  the  case  of  orthographic  projection,  one  constraint  on  and  was  that  they  each  have 
unit  magnitude,  as  required  by  equation  (10).  In  the  above  paraperspective  version  of  those 
constraints,  we  simply  require  that  their  magnitudes  be  in  a  certain  ratio. 


There  is  also  a  constraint  on  the  angle  relationship  of  and  From  the  definition  of 
and  Up  we  have 


CX^Jf 


(31) 
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The  problem  with  this  constraint  is  that,  again,  is  unknown.  We  could  choose  to  use  either 
value  from  equation  (29)  for  l/zj,  since  theoretically  they  should  be  equal,  but  in  order  to 
avoid  relying  on  some  values  more  than  others,  we  use  the  average  of  the  two  quantities.  We 
choose  the  arithmetic  mean  over  the  geometric  mean  or  some  other  measure  in  order  to  keep 
the  constraints  linear  mQ  =  A^A.  Thus  our  second  constraint  becomes 

=  (32, 

^  ^  2(1  + exp  2(\  +  cyp 

This  is  the  paraperspective  version  of  the  orthographic  constraint  given  by  equation  (11), 
which  required  that  the  dot  product  of  and  ii^  be  zero. 

Equations  (30)  and  (32)  are  homogeneous  constraints,  which  could  be  trivially  satisfied  by 
the  solution  M  =  0.  To  avoid  this  solution,  we  impose  the  additional  constraint  that 

h,|  =  1.  (33) 

This  does  not  effect  the  final  solution  except  by  a  scaling  factor. 

Equations  (30),  (32),  and  (33)  are  the  paraperspective  version  of  the  metric  constraints,  and 
we  compute  the  3  x  3  matrix  A  such  that  M  =  MA  best  satisfies  the  metric  constraints  in  the 
least  sum-of-squares  error  sense.  This  is  a  simple  problem  because  we  have  been  careful  to 
ensure  that  they  are  linear  constraints  in  the  6  unique  elements  of  the  symmetric  3x3  matrix 
Q  =  A.  We  use  the  metric  constraints  to  compute  Q,  compute  its  Jacobi  Transformation 
Q  ~  lmJ,  where  A  is  the  diagonal  eigenvalue  matrix,  and  as  long  as  2  is  positive  definite, 
A  is 


A  =  (LA^)^.  (34) 

4.4.  Shape  and  Motion  Recovery 

Once  the  matrix  A  has  been  determined,  the  shape  matrix  S  and  the  motion  matrix  M  can 
easily  be  computed  fitom  equation  (27).  For  each  frame  f,  however,  there  is  a  more  complex 
relationship  between  the  actual  translation  and  rotation  vectors  and  the  rows  of  the  matrix 
M.  From  equation  (19)  we  can  see  that 

V=V“/+cxyt/  ]f=zjaf+cy/if.  (35) 

We  use  this  to  determine  the  following  equations: 

(tf^f+cx^f)  X  (zja^+cyjkf)  =kf 

III  =  I  */*»/+ cJfAl  =  1  (36) 

\]i  =  \z^f+cyjk^  =  l 

Again,  we  do  not  know  a  value  for  zp  but  using  the  relations  specified  in  equation  (29)  and 


the  additional  knowledge  that  fij  »  1 ,  equation  (36)  can  be  reduced  to 
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(i/xo/)  •&/=  1 

m^t/=  -cXf  (37) 

nfkf=-cyf 


where  m#  =  ]\  +  cxl^  and  if  -  Jl  +  cyi^.  These  constraints  are  put  in  the  form 


G^f  »  Hf, 


(38) 


where 


■(myrxi/)' 

-  1  - 

“/ 

"/  = 

-CXf 

-  “/  - 

We  confute  fe/  simply  as 


kf  =  g;*//, 

and  then  compute 


(39) 


(40) 


If^ifxif  ]f=kfxmf.  (41) 

Because  and  may  not  have  exactly  satisfied  the  metric  constraints,  there  is  no  guaran¬ 
tee  that  the  1/  and  j/  given  by  this  equation  will  be  orthonormal.  Therefore  we  actually  com¬ 
pute  the  orthonormal  1/  and  j/  which  are  closest  to  the  values  given  by  equation  (41).  Due  to 
the  arbitrary  world  coordinate  orientation,  to  obtain  a  unique  solution  we  then  rotate  the 
computed  shape  and  motion  to  align  the  world  axes  with  the  first  firame’s  camera  axes,  so 
that  ii  =  [i  0  and  ji  =  [o  i 

All  that  remain  to  be  computed  are  the  translations  for  each  frame.  We  calculate  the  depth 
from  either  part  or  some  combination  of  the  parts  of  equation  (29).  Since  we  already  know 
cxp  cyp  1/  and  j/,  we  calculate  t/  using  equations  (17)  and  (18). 


4.5.  Solution  Ambiguity  Removal 


In  order  to  solve  for  A  from  Q,  in  equation  (34)  we  took  the  square  root  of  the  diagonal 


eigenvalue  noatrix  A.  It  is  clear,  however,  that  any  matrix  of  the  form  A 


±10  0 
0  ±1  0 
0  0  ±1 


would 


produce  the  same  matrix  Q  when  multiplied  by  its  transpose.  This  sign  ambiguity  in  the  first 
two  columns  of  A  was  removed  when  we  aligned  the  world  coordinate  axes  with  the  first 
frame’s  camera  axes,  at  the  end  of  Section  4.4.  However,  the  ambiguity  in  the  third  column 
of  A  remains.  This  affects  the  final  shape  simply  by  negating  the  z-component,  and  effects 
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the  rotation  by  negating  the  third  column  of  the  matrix  M. 

Geonoetrically,  this  ambiguity  arises  because  paraperspecdve  projection  does  not  account 
for  any  different  perspective  deformation  within  the  object  itself.  Thus  it  is  not  possible  to 
distinguish  the  “tont”  of  the  object  from  the  “back”  of  the  object,  as  can  be  seen  from  Fig¬ 
ure  5(a).  However,  in  real  scenarios  there  will  be  additional  queues  due  to  occlusion  infor¬ 
mation  as  in  Figure  5(b)  and  (c),  or  due  to  perspective  distortion  of  the  object,  as  in  Figure 
5(d),  if  the  object  is  not  too  distant  from  the  canwra.  Simple  methods  based  on  either  of 
these  phenomena  should  be  sufficient  to  determine  which  of  the  two  solutions  given  by  the 
paraperspecdve  factorizadon  method  is  consistent  with  the  image  data. 


<*>  (b)  (c)  (d) 

Figure  5  Ambiguity  of  Soiution 

(a)  Sequence  of  images  with  two  valid  motion  and  slugw  interpretations. 

(c)  Ambiguity  removed  due  to  occlusion  infmmation  in  the  image  sequence. 

(d)  Ambiguity  removed  due  to  perq)ective  distortion  of  the  object  in  the  images 

4.6.  Paraperspective  Projection  as  an  Approximation  to  Perspective 
Projection 

In  Secdon  4.1.,  we  defined  paraperspecdve  projecdcm  geometrically.  We  can  arrive  at  the 
same  equadons  mathematically  as  a  first-order  approximadon  to  the  perspecdve  projecdon 
equadons.  The  perspecdve  projecdon  equadons  are 


where 


Vp 


^/p  “  *'/■  ^*p“9 

For  simplicity  we  assume  unit  focal  length,  /  =  i . 

Taking  the  Taylor  series  expansion  of  these  equations  about  the  point 


yields 


Vp-V=-V*V 


(43) 


(44) 


_  /  ^  P  T  f  '  P  “  P 

V  ^ 

‘f  ' 


(2fp-Zf)+2^ 


■) 


-V’(®P  9  r,  _,)2 
^fP  9 


(Zfp-Zf)  +2^ - 


4 


(45) 


Ignoring  all  but  the  first  term  of  the  Taylor  series  yields  the  equations  for  scaled  orthogra¬ 
phic  projection  (Sec  Appendix  I.)  However,  instead  of  arbitrarily  stopping  at  the  first  term, 
we  can  eliminate  higher  order  terms  based  on  the  approximation  that  =  0.  Using  this 

approximation,  and  combining  equations  (43)  and  (44)  to  determine  that  z/p-Zf  =1^  5^, 
gives  the  following  equations 


T  -  */  9  V  <  9  . 

9p - T  2  V 

^  V 

_i/  <»p-V  h  <-V) 


(46) 


>/P  = 


(k/V 


Simply  factoring  out  the  l/z,  and  expanding  the  dot-products  i^p-i/i  and  j^-  (s^-ty) 
gives 


^fp = riv-  >p-"  V  -  V> ) 

^  I  •  t 

y/p = *P + V  -  V> ) 

This  is  equivalent  to  the  paraperspective  projection  equations  given  by  equation  (15). 

The  approximation  that  =  0  preserves  a  portion  of  the  second  term  of  the  Taylor 

series  expansion,  the  term  which  is  of  order  (js^tyj)  /z^,  while  ignoring  the  remaining  por¬ 
tion  of  the  second  term,  which  is  of  order  |Sp| and  all  higher  order  terms.  Qearly  if  the 
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translation  that  the  object  undergoes  is  also  small,  then  there  is  little  justification  for  preserv¬ 
ing  this  portion  of  the  second  term  and  not  the  other.  In  such  cases,  the  entire  second  term 
can  be  safely  ignored,  leaving  only  the  equations  for  scaled  orthographic  projection. 

One  will  notice  that  we  did  not  explicitly  set  the  world  origin  at  the  object’s  center-of-mass 
as  we  did  in  Section  4.1.  However,  the  assumption  that  |Sp|^/r/  =  0  will  be  most  accurate 
when  the  magnitudes  of  the  are  smallest.  Since  the  vectors  represent  the  vectors  from 
the  wt^d  origin  to  the  object  points,  their  magnitudes  will  be  snrallest  when  the  world  ori¬ 
gin  is  located  at  the  object’s  center-of-mass. 
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5.  Experiments 


5.1.  Parameters 

To  test  our  paraperspective  factorization  toethod,  we  created  synthetic  point  sequences  using 
a  perspective  projection  model  of  objects  undergoing  rotation  and  translation.  We  then  per¬ 
turbed  the  coordinates  of  each  point  by  adding  gaussian  noise.  We  used  three  different 
object  shapes;  a  box  object  whose  edges  were  of  unit  length,  a  sphere  with  unit  diameter, 
and  a  random  point  cloud  within  a  unit  cube,  each  containing  approximately  60  feature 
points.  All  of  the  test  runs  consisted  of  60  image  frames  of  the  object  rotating  through  a  total 
of  30  degrees  each  of  roll,  pitch,  and  yaw. 

Our  goal  was  not  only  to  ascertain  the  performance  of  the  paraperspective  factorization 
method,  but  to  determine  its  performance  in  comparison  to  other  methods.  We  used  five  dif¬ 
ferent  methods  to  recover  the  shape  and  motion  of  the  object  and  compared  the  accuracy  of 
the  results.  The  first  method  was  the  orthographic  factorization  metho<L  The  second  was  the 
scaled  orthographic  factorization  method,  which  is  described  in  Appendix  I.  The  third 
method  was  our  paraperspective  factorization  method.  The  fourth  method  was  a  full  per¬ 
spective  method,  which  iteratively  solves  the  perspective  projection  equations,  as  described 
in  Appendix  H.  For  this  method,  we  used  the  paraperspective  factorization  method  to  deter¬ 
mine  an  initial  value.  Because  such  methods  are  generally  very  sensitive  to  initial  condi¬ 
tions,  our  fifth  method  used  the  same  iterative  solution  technique,  but  we  used  the  true  shape 
and  motion  as  initial  values.  This  method  is  unfortunately  not  an  option  in  real  systems,  but 
indicates  what  is  essentially  an  upper  bound  on  the  accuracy  achievable  using  a  least  sum- 
of-squares  difference  formulation  of  the  full  perspective  projection  model,  without  making 
any  assumptions  about  the  motion  or  object  shape. 

In  particular,  we  were  interested  in  how  the  methods  compare  as  a  function  of  the  noise  level 
in  the  image  and  the  distance  of  the  object  from  the  camera.  Thus  the  depth,  representing  the 
distance  from  the  camera’s  focal  point  to  the  front  of  the  object  in  the  first  fi-ame,  was  varied 
from  3  to  60  times  the  object  size.  In  generating  our  synthetic  images,  for  each  depth  we 
chose  the  largest  focal  len^  which  would  keep  the  object  in  the  field  of  view.  We  varied  the 
standard  deviation  of  the  gaussian  noise  from  0  to  4  pixels  (of  a  512x512  pixel  image)  stan¬ 
dard  deviation.  For  each  combination  of  object,  depA,  and  noise,  we  performed  three  tests, 
using  different  random  noise  each  time,  and  averaged  the  resulting  errors. 

5.2.  Error  Measurement 

We  present  here  the  total  shape  error,  rotation  error,  X-Y  offset  error,  and  Z  offset  (depth) 
error.  The  term  “offset”  refers  to  translations  along  the  camera  components;  the  X  offset  is 
\f  \f,  the  Y  offset  is  1/  -  3/.  and  the  Z  offset  is  t/-  kf.  These  are  the  components  of  t/  actually 
recovered  by  the  methods,  and  although  a  simple  transformation  could  recover  the  three  ele¬ 
ments  of  \f  in  world  coordinates,  the  resulting  translation  errors  would  depend  on  the  accu¬ 
racy  of  the  recovered  rotation.  Note  that  the  shape  and  translation  can  only  be  determined  up 
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to  a  scaling  factor  (it’s  not  possible  to  tell  whether  the  image  sequence  is  of  a  house  fifty 
meters  away  or  of  a  1/10  scale  model  of  a  house  five  meters  away).  To  compute  the  shape 
error,  we  find  a  scale  which  minimizes  the  root-mean-squaie  (RMS)  error  between  the  true 
and  computed  shape,  and  then  return  this  error.  We  use  the  same  method  for  the  X-Y  offset, 
and  for  the  Z  offset  This  does  not  correspond  exactly  to  a  valid  solution,  because  in  reality 
there  is  only  a  single  scaling  factor  for  all  of  the  shape  and  offset  values,  but  this  avoids 
cases  in  wMch,  for  example,  a  large  shape  error  influences  the  scaling  factor,  which  then 
causes  the  offset  error  to  ht  reported  as  large.  The  rotation  error  is  computed  as  the  RMS  of 
the  size  in  radians  of  the  angle  by  which  a  computed  camera  frame  must  be  rotated  about 
some  axis  to  produce  the  tme  camera  frame. 

5.3.  Results 

In  the  experiments  in  which  the  object  was  centered  in  the  image  and  there  was  no  depth 
translation,  we  found  that  the  orthographic  factorization  method  performed  well,  and  the 
paraperspective  factorization  method  provided  no  significant  improvement.  However,  in 
image  sequences  in  which  there  was  depth  translation  or  the  object  was  not  centered  in  the 
image,  the  paraperspective  method  performed  significantly  better  than  the  orthographic  fac¬ 
torization  method.  The  results  of  these  experiments  are  shown  on  the  following  pages.  The 
average  error  results  were  very  similar  for  all  of  the  objects,  so  our  graphs  show  the  average 
over  all  runs  of  all  objects.  Since  there  were  3  objects  and  we  performed  3  runs  per  object 
with  different  random  noise,  each  data  point  represents  the  average  of  9  runs. 

Figure  6  shows  how  the  various  methods  performed.  In  these  image  sequences,  the  object 
moved  across  the  screen  one  unit  horizontally  and  one  unit  vertically.  The  object  also  moved 
away  from  the  camera  by  a  total  of  one  half  the  object’s  initial  distance  from  the  camera. 
Thus  in  a  test  case  in  which  the  object’s  depth  in  the  first  frame  was  3.0,  its  depth  in  the  last 
frame  was  4.5.  These  experiments  were  done  using  a  noise  standard  deviation  of  2  pixels, 
which  we  consider  a  rather  high  noise  level.  At  low  depths,  perspective  distortion  is  a  signif¬ 
icant  source  of  error  in  the  computed  results.  Interestingly,  our  experiments  show  that  for 
objects  farther  from  the  camera  than  7  tiroes  the  object  size,  refining  the  paraperspective 
solution  using  the  perspective  iteration  technique  improves  the  rotation  and  translation  very 
little.  However,  even  at  depths  beyond  30  times  the  object  size,  the  perspective  refinement 
method  significantly  improves  the  shape. 

The  behavior  of  the  paraperspective  factorization  method  over  a  range  of  noise  levels  is 
shown  in  Figure  7.  We  see  that  once  the  object  is  far  enough  firom  the  camera  that  perspec¬ 
tive  effects  are  minor,  the  error  in  the  computed  solution  is  nearly  proportional  to  the 
amount  of  noise  in  the  input  For  some  reason  which  we  are  currently  unable  to  explain,  the 
rotation  oror  seems  to  vary  widely  even  at  a  given  noise  level. 

We  implemented  the  methods  in  C  and  performed  the  experiments  on  a  Sun  4/65.  Solving 
the  system  of  60  frames  and  60  points  required  about  20-24  seconds  for  each  of  the  three 
factorization  methods,  with  much  of  this  time  spent  computing  the  singular  value  decompo¬ 
sition  of  the  measureoMnt  matrix.  The  perspective  iteration  method,  however,  required 
about  250  seconds  to  solve  the  same  system. 


X  and  Y  Offset  Enor  X  lOc-3  Rotation  Error  x  lOe-3 
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F^ure  6  Comparison  of  Methods  for  a  Topical  Case 
These  figures  diow  the  error  in  the  computed  solution  for  the  five  algorithms,  as  a  function  of  the 
object’s  distance  fiom  the  camera.  The  image  sequences  contained  translation  across  the  image 
and  away  from  the  camera,  and  contained  added  gaussian  noise  of  2  pixels  standard  deviation. 


X  and  Y  Offset  Error  x  lOc-6  Rotation  Error  x  lOe-3 
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Future  7  Parapo-spective  Shape  md  Motion  Recovery  by  Noise  Level 
These  figures  diow  the  error  in  the  solution  computed  by  die  paraperspective  factorization 
method  for  a  range  of  different  noise  levels,  as  a  function  of  the  object’s  distance  from  the  cam¬ 
era.  These  sequences  contained  translation  across  the  image  and  away  from  the  camera. 
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6.  Conclusions 

The  principle  that  the  measurement  matrix  has  rank  3,  as  put  forth  by  Tomasi  and  Kanade  in 
[7],  depended  on  the  use  of  an  orthographic  projecdon  model.  We  have  shown  in  this  paper 
that  this  important  result  also  holds  for  ^e  case  of  paraperspective  projection,  which  closely 
approximates  perspective  projection,  and  have  devised  a  paraperspective  factorization 
method  based  on  this  projection  model. 

In  general  image  sequences  in  which  the  object  being  viewed  translates  significantly  toward 
or  away  from  the  camera  or  across  the  camera’s  field  of  view,  the  paraperspective  method 
performs  significantly  better  than  the  method  based  on  orthographic  projection.  In  image 
sequences  in  which  the  orthographic  assumption  is  valid,  where  the  object  is  centered  in  the 
image  and  not  translating  along  the  camera’s  optical  axis,  the  orthographic  factorization 
method  often  produces  better  results  than  the  paraperspective  or  even  the  full  perspective 
methods.  This  is  because  the  orthographic  method  is  in  effect  using  known  information 
about  the  motion  in  the  sequence,  adding  the  constraint  that  the  object  does  not  move  toward 
or  away  firam  the  camera.  In  images  sequences  in  which  the  object  is  close  to  the  camera,  the 
paraperspective  factorization  method  still  provides  accurate  motion  results,  and  provides  a 
good  approximation  of  the  object  shape;  this  solution  can  be  further  refined  using  an  itera¬ 
tive  perspective  method. 

The  paraperspective  factorization  method  computes  the  distance  from  the  camera  to  the 
object,  which  enables  its  use  in  a  wider  range  of  scenarios.  The  method  performs  well  in  a 
wide  range  of  motion  scenarios,  is  efficient,  and  is  robust  with  respect  to  noise. 

In  real  image  sequences,  feature  points  will  often  become  occluded  and  new  feature  points 
will  appear.  This  problem  of  occlusion  has  been  handled  for  the  orthographic  factorization 
method  [8],  but  has  not  been  addressed  here  for  the  paraperspective  factorization  method. 
We  will  similarly  extend  the  paraperspective  factorization  method  to  accommodate  occlu¬ 
sion. 

So  far  we  have  concentrated  on  testing  the  method  on  synthetic  image  sequences  in  order  to 
observe  its  performance  across  a  very  wide  range  of  situations  and  noise  levels,  and  have 
performed  only  preliminary  testing  on  real  image  sequences.  In  order  to  verify  the  method’s 
practicality,  we  plan  to  perform  more  extensive  tests  on  real  image  sequences. 
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Appendix  I.  Scaled  Orthographic  Factorization  Method 

Scaled  orthographic  projection,  also  known  as  “weak  perspective”  [3],  is  a  closer  approxi¬ 
mation  to  perspective  projection  than  orthographic  projection,  yet  not  as  accurate  as  parap- 
erspective  projection.  This  method  can  be  used  when  the  object  remains  centered  in  the 
image,  or  when  the  distance  m  the  object  is  large  relative  to  the  size  of  the  object. 

1.1.  Scaled  Orthographic  Projection  Equations 

In  the  scaled  orthographic  projection  model,  object  points  are  projected  along  an  axis  paral¬ 
lel  to  the  camera’s  optical  axis,  onto  a  plane  perpendicular  to  the  optical  axis  and  passing 
through  the  object’s  center-of-mass  c.  This  image  is  then  projected  onto  the  image  plane 
using  perspective  projection,  but  because  the  points  are  on  a  plane  parallel  to  the  image 
plane,  this  is  equivalent  to  simply  scaling  the  points  by  the  distance  along  the  optical  axis 
from  the  camera  to  the  object’s  center-of-mass. 

Image 


Scaled  Orthographic  Projection  in  two  dimensions 
Dotted  lines  indict  trae  perspective  projection 

This  changes  the  orthographic  projection  equations  by  adding  a  scaling  factor  dependent  on 
the  distance  to  the  object’s  center-of-mass,  which  is  given  by 

(c-y 

Thus  the  scaled  orthographic  projection  equations  are 


(48) 
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^/P  =  7//- 

y/p=7^^f  <»p-v» 


(49) 


To  simplify  the  equations  we  assume  unit  focal  length,  /  =  i .  Because  the  world  origin  is 
arbitrary,  we  set  the  world  origin  at  the  object’s  center-of-mass,  so  that  c  =  0,  and  rewrite 
the  above  equations  as 


where 


=  mpip  +  cxy 

yfp=nfBp^cyp 

(50) 

It 

-V  »'/ 

(51) 

V‘  */ 

CXf  =  -- — i 

*/■■•/ 

=  — 

(52) 

(53) 

1.2.  Decomposition 

Because  equation  (SO)  is  identical  to  equation  (2),  the  measurement  matrix  w  can  still  be 
written  as  W'  «  ATS  +  r  just  as  in  orthographic  case.  We  can  still  compute  cx^  and  cy.  immedi¬ 
ately  from  the  image  data  using  equation  (24),  and  use  singular  value  decomposition  to  fac¬ 
tor  the  registered  measurement  matrix  W  into  the  product  of  ik  and 

IJ.  Normalization 

Again,  the  decomposition  is  not  unique  and  we  must  determine  the  3  x  3  matrix  A  which 
produces  the  actual  motion  matrix  M  =  (kA  and  the  shape  matrix  S  =  A~^S.  We  see  from 
equation  (53)  that 


1“/  =  \  I"/  =  (54) 

We  do  not  know  the  value  of  the  depth  Zp  so  we  cannot  impose  individual  constraints  on 
and  as  we  did  in  the  orthographic  case.  Instead,  we  combine  the  two  equations  as  we  did 
in  the  paraperspective  case,  to  impose  the  following  constraint; 

1“/  =  !“/•  (55) 

Because  and  are  just  scalar  multiples  of  and  jp  we  can  still  use  the  constraint  that 


in^n^  =  0. 


(56) 
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As  in  the  paraperspective  case,  equations  (SS)  and  (56)  are  homogeneous  constraints,  which 
could  be  trivially  satisfied  by  the  solution  M  =  0,  so  to  avoid  this  solution  we  add  the  some¬ 
what  arbitrary  consumnt  that 

imil  =  1.  (57) 

Equations  (55),  (56),  and  (57)  are  the  scaled  orthographic  version  of  the  metric  constraints. 
We  can  compute  the  3  x  3  matrix  A  which  best  satisfies  them  very  easily,  because  the  con¬ 
straints  are  linear  in  the  6  unique  elements  of  the  symmetric  3x3  matrix  Q  =  aJa. 

1.4.  Shape  and  Motion  Recovery 

Once  the  matrix  A  has  been  found,  the  shape  can  easily  be  computed  from  equation  (9).  We 
compute  the  motion  parameters  by 

Unlike  the  orthographic  case,  we  can  now  compute  the  depth  to  the  object  center-of-mass  in 
each  firame,  using  either  expression  in  equation  (54),  or  by  using  the  geometric  mean  of  the 
two  solutions. 


(59) 


Appendix  U.  Perspective  Method 


This  section  presents  the  iterative  method  used  to  recover  the  shape  and  motion  using  a  per¬ 
spective  projection  model.  Although  our  algorithm  was  developed  independently  and  han¬ 
dles  the  full  three  dimensional  case,  this  method  is  quite  similar  to  a  two  dimensional 
algorithm  developed  by  Taylor,  Kriegman,  and  Anandan  as  reported  in  [6]. 

n.l.  Perspective  Projection  Equations 

In  the  perspective  projection  model,  sometimes  referred  to  as  the  pinhole  camera  model, 
each  point  is  projected  onto  the  image  plane  directly  towards  the  focal  point  of  the  camera. 
An  object  point’s  image  coordinates  are  determined  by  the  position  at  which  the  line  con¬ 
necting  the  object  point  with  the  camera’s  focal  point  intersects  the  image  plane.  This  is 
illustrated  in  Figure  9. 


Image 


Simple  geometry  using  similar  triangles,  gives  us 
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We  now  assume  unit  focal  length  and  rewrite  this  as 


fP 

^fp  VSp  +  «/ 

where 

«/=-*/■*/  ''>’/  = -vv 

n.2.  Iterative  Minimization  Method 


“/  = 


V 


(61) 


(62) 


For  this  method,  we  take  advantage  of  the  fact  that  y,  and  can  be  written  as  functions 
of  only  three  independent  rotational  parameters. 


cosa^os^^  (cosaySinP^sinT^-sinayCOST^)  (cosaySinpyCosT^H- sinUySinT^) 
sinOyCOsPy  (sinaySinpySm-iy+cosayCOSTy)  (sinaySinPyCOSTiy- cosOySinTy) 
-sinPy  cosPySinYy  cosPyCOS')^ 


(63) 


Thus  there  are  six  motion  parameters,  cxy,  cyy,  czp  Uy,  Py,  and  for  each  firame,  and  three 
sha^  parameters,  =  [jp^j  f  ^  for  each  point.  Equation  (61)  gives  us  an  overcon¬ 
strained  set  of  2FP  equations  m  these  6F  +  3P  variables.  Because  the  perspective  projection 
equations  are  non-linear,  there  is  no  apparent  way  to  combine  the  equations  for  all  points 
and  frames  into  a  single  matrix  equation,  as  we  did  in  the  orthographic  and  paraperspective 
case,  or  to  separate  the  shape  from  the  motion,  or  to  compute  the  translational  components 
directly  from  the  image  measurements.  Instead,  we  formulate  the  problem  as  a  general  least- 
squares  problem  in  which  we  want  to  determine  the  values  of  the  6F+  3P  variables  which 
minimize 


(64) 


We  could  in  theory  apply  any  one  of  a  number  of  non-linear  equation  solution  techniques  to 
this  problem.  Such  methods  begin  with  a  set  of  initial  starting  values,  and  iteratively  refine 
those  values  to  reduce  the  error.  Because  we  know  the  form  of  the  equations,  we  can  use 
derivative  information  to  guide  our  numerical  search.  However,  general  non-linear  least 
square  techniques  would  not  take  full  advantage  of  the  structure  of  our  equations.  For  every 
frame  f  and  point  p,  we  would  have  to  compute  (dxyp  /  (dv)  and  (dyy^)  /  (9v)  for  each  of  the 
6F  +  3P  variables  v.  However,  in  our  equations  most  of  these  derivatives  are  zero.  In  fact. 
For  a  given  frame  f  and  point  p,  these  partial  derivatives  are  non-zero  only  for  eight  of  the 
6F  +  3P  variables. 


Our  solution  technique  takes  advantage  of  the  particular  structure  of  the  equations  by  sepa¬ 
rately  refining  the  shape  and  motion  parameters.  We  hold  the  shape  constant  and  solve  for 
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the  motion  parameters  which  minimize  the  error.  We  then  hold  the  motion  constant,  and 
solve  for  the  shape  parameters  which  minimize  the  error.  We  repeat  this  process  until  an  iter¬ 
ation  produces  no  significant  reduction  in  the  total  sum  of  squares  error.  While  holding  the 
shape  constant,  the  minimization  can  be  done  independently  for  each  firame,  which  involves 
minimizing  the  error  in  a  system  of  6  variables  in  P  equations.  Likewise  wliile  holding  the 
motion  constant,  we  can  solve  for  the  shape  separately  for  each  point  by  solving  a  system  of 
2F  equations  in  3  variables.  This  not  only  reduces  the  problem  to  manageable  complexity, 
but  as  pointed  out  in  [6],  it  lends  itself  well  to  parallel  implementation. 

We  perform  the  individual  minimizations,  fitting  six  motion  variables  to  P  equations  or  fit¬ 
ting  three  shape  variables  to  IF  equations,  using  the  Levenberg-Marquardt  method  [S],  a 
method  which  uses  steepest  descent  when  far  from  the  minimum  and  varies  continuously 
towards  the  inverse-Hessian  method  as  the  minimum  is  approached. 

In  fact  we  do  not  need  to  vary  all  6F+ZP  variables.  Recall  that  the  solution  is  only  deter¬ 
mined  up  to  a  scaling  factor,  the  world  origin  is  arbitrary,  and  the  world  coordinate  orienta¬ 
tion  is  arbitrary.  These  factors  reduce  the  number  of  variable  assignments  representing 
different  real  solutions.  We  could  choose  to  arbitrarily  keep  the  first  frame’s  rotation  fix  at 
zero  roll,  pitch,  and  yaw,  and  similarly  fix  some  shape  or  translation  parameters.  However, 
we  found  in  our  experiments  that  the  algorithm  converged  significantly  faster  if  all  of  the 
shape  and  motion  parameters  were  allowed  to  vary.  Once  the  algorithm  converges  to  a  solu¬ 
tion,  we  adjust  the  final  shape  and  translation  to  place  the  origin  at  the  object’s  center-of- 
mass,  scale  the  solution  so  that  the  depth  in  the  first  frame  is  1,  and  rotate  the  solution  so  that 
*1  =  [l  0  ^  and  ji  =  [o  1  ^^,orequivalently  in  this  formulation,  so  that  Oj  =  Pj  =  Fj  =  0. 

One  problem  with  this  sort  of  iteration  method  is  that  its  final  result  can  be  highly  dependent 
on  its  initial  values.  Taylor,  Kriegman,  and  Anandan  [6]  require  some  basic  odometry  mea¬ 
surements  as  might  be  produced  by  a  navigation  system  to  use  as  initial  values  for  their 
motion  parameters,  and  use  the  2D  shape  of  the  object  in  the  first  image  firame,  assuming 
constant  depth,  as  their  initial  shape.  To  avoid  the  requirement  for  odometry  measurements, 
which  way  not  be  available  in  many  situations,  we  use  the  paraperspective  factorization 
method  to  supply  initial  values  to  the  perspective  iteration  method. 


